#https://hvplot.holoviz.org/user_guide/Gridded_Data.html
# Curvilinear http://holoviews.org/user_guide/Gridded_Datasets.html
import os
import xarray as xr
import numpy as np
import pandas as pd

from utils.misc_utils import is2_interp2d, restrictRegionally
from utils.read_data_utils import read_is2_data, read_book_data

# Plotting
import hvplot.xarray
import cartopy.crs as ccrs
#import cartopy.feature as cfeature
#from matplotlib.axes import Axes
#from cartopy.mpl.geoaxes import GeoAxes
#from textwrap import wrap
#import matplotlib as mpl
#import matplotlib.pyplot as plt

#GeoAxes._pcolormesh_patched = Axes.pcolormesh
#%config InlineBackend.figure_format = 'retina'
#plt.rcParams['axes.facecolor']='white' 
#mpl.rcParams['figure.dpi'] = 200
# Set winter months 
winter18_19 = pd.date_range(start="2018-11-01", end="2019-04-01", freq="MS")
winter19_20 = pd.date_range(start="2019-09-01", end="2020-04-01", freq="MS")
winter20_21 = pd.date_range(start="2020-09-01", end="2021-04-01", freq="MS") 
winter_months = winter18_19.append(winter19_20).append(winter20_21)

ds_uninterpolated = read_is2_data() # Read in data

# Read book data. This contains the CDR data, which we'll use for interpolating 
book_ds = read_book_data()
book_ds = book_ds.sel(time=winter_months) # Get winter months

# Interpolate
print("Interpolating ICESat-2 data...")
cdr_da = book_ds["cdr_seaice_conc_monthly"] # Get CDR data
ds_interpolated = is2_interp2d(ds_uninterpolated, cdr_da, method='nearest', interp_var='all')
print("Complete!")
    
# Combine datasets 
ds_interpolated = xr.merge([ds_interpolated, book_ds])
ds_interpolated = ds_interpolated.drop_vars("projection")
ds_uninterpolated = ds_uninterpolated.drop_vars("projection")
Interpolating ICESat-2 data...
Complete!
---------------------------------------------------------------------------
MergeError                                Traceback (most recent call last)
/var/folders/8v/cr06mz0n3bjd5jr12_9v6n200000gn/T/ipykernel_81870/369422530.py in <module>
     18 
     19 # Combine datasets
---> 20 ds_interpolated = xr.merge([ds_interpolated, book_ds])
     21 ds_interpolated = ds_interpolated.drop_vars("projection")
     22 ds_uninterpolated = ds_uninterpolated.drop_vars("projection")

~/opt/anaconda3/envs/icesat2_book/lib/python3.8/site-packages/xarray/core/merge.py in merge(objects, compat, join, fill_value, combine_attrs)
    898         dict_like_objects.append(obj)
    899 
--> 900     merge_result = merge_core(
    901         dict_like_objects,
    902         compat,

~/opt/anaconda3/envs/icesat2_book/lib/python3.8/site-packages/xarray/core/merge.py in merge_core(objects, compat, join, combine_attrs, priority_arg, explicit_coords, indexes, fill_value)
    633 
    634     prioritized = _get_priority_vars_and_indexes(aligned, priority_arg, compat=compat)
--> 635     variables, out_indexes = merge_collected(
    636         collected, prioritized, compat=compat, combine_attrs=combine_attrs
    637     )

~/opt/anaconda3/envs/icesat2_book/lib/python3.8/site-packages/xarray/core/merge.py in merge_collected(grouped, prioritized, compat, combine_attrs)
    238                 variables = [variable for variable, _ in elements_list]
    239                 try:
--> 240                     merged_vars[name] = unique_variable(name, variables, compat)
    241                 except MergeError:
    242                     if compat != "minimal":

~/opt/anaconda3/envs/icesat2_book/lib/python3.8/site-packages/xarray/core/merge.py in unique_variable(name, variables, compat, equals)
    147 
    148     if not equals:
--> 149         raise MergeError(
    150             f"conflicting values for variable {name!r} on objects to be combined. "
    151             "You can skip this check by specifying compat='override'."

MergeError: conflicting values for variable 'ice_thickness_smoothed' on objects to be combined. You can skip this check by specifying compat='override'.
#https://hvplot.holoviz.org/user_guide/Geographic_Data.html
date="2021-04"
thickness_interp = ds_interpolated["ice_thickness"]
thickness_uninterp = ds_uninterpolated["ice_thickness"]

# PLOT!! 
thickness_uninterp.sel(time=date)[0].hvplot.quadmesh(y="latitude",x="longitude", 
                                                colorbar=False, clim=(0,4), clabel="sea ice thickness (m)",
                                                projection=ccrs.NorthPolarStereo(central_longitude=-45), cmap="viridis", features=["coastline"], 
                                                title="uninterpolated sea ice thickness "+pd.to_datetime(date).strftime("%Y-%m"), 
                                                project=True, geo=True, width=450, height=300, ylim=(60,90)) + \
thickness_interp.sel(time=date)[0].hvplot.quadmesh(y="latitude",x="longitude", 
                                              colorbar=True, clim=(0,4), clabel="sea ice thickness (m)",
                                              projection=ccrs.NorthPolarStereo(central_longitude=-45), cmap="viridis", features=["coastline"], 
                                              title="interpolated sea ice thickness "+pd.to_datetime(date).strftime("%Y-%m"), 
                                              project=True, geo=True, width=500, height=300, ylim=(60,90))
thickness_interp = ds_interpolated["ice_thickness"]
date="2021-04"
thickness_interp.sel(time=date)[0].hvplot.quadmesh(y="latitude",x="longitude", 
                                                colorbar=False, clim=(0,4), clabel="sea ice thickness (m)",
                                                projection=ccrs.NorthPolarStereo(central_longitude=-45), cmap="viridis", features=["coastline"], 
                                                title="uninterpolated sea ice thickness "+pd.to_datetime(date).strftime("%Y-%m"), 
                                                project=True, geo=True, width=450, height=300, ylim=(60,90))